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We study the stability of three-dimensional incompressible Weyl semimetals in the presence of 
random quenched charge impurities. Combining numerical analysis and scaling theory we show 
that in the presence of sufficiently weak randomness (i) Weyl semimetal remains stable, while (ii) 
double-Weyl semimetal gives rise to compressible diffusive metal where the mean density of states 
at zero energy is finite. At stronger disorder, Weyl semimetal undergoes a quantum phase transition 
and enter into a metallic phase. Mean density of states at zero energy serves as the order parameter 
and displays single-parameter scaling across such disorder driven quantum phase transition. We 
numerically determine various exponents at the critical point, which appear to be insensitive to the 
number of Weyl pairs. We also extract the extent of the quantum critical regime in disordered Weyl 
semimetal and the phase diagram of dirty double Weyl semimetal at finite energies. 
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Introduction : Over the span of last few years the hori¬ 
zon of topological phases of matter has been extended 
beyond the gapped states [Hi and now includes various 
nodal (gapless) systems as well [3; 0] . Three dimensional 
Weyl semimetal (WSM) is the prime example of such 
non-insulating systems, which is constituted by so called 
Weyl nodes that act as source (monopole) and sink (anti¬ 
monopoles) for Berry flux in the reciprocal space, thus 
always appear in pairs [5]. As the hallmark signature of a 
topologically nontrivial phase, WSM accommodates gap¬ 
less (chiral) surface states that can give rise to peculiar 
electro-magnetic responses, such as, anomalous Hall and 
chiral-magnetic effects [B]. While the gapped topologi¬ 
cal phases are expected to be robust against sufficiently 
weak randomness, stability of their gapless counterpart 
against disorder demands careful investigation and con¬ 
stitutes the central theme of this Rapid Communication. 

Recent time has witnessed the discovery of WSMs in 
a number of noncentrosymmetric and magnetic semi¬ 
conductors IMU. Various other proposals for WSMs, 
for example, include anti-ferromagnetically m or spin- 
ice j!8j ordered pyrochlore iridates, multilayer configura¬ 
tion of topological and regular insulators mum, mag¬ 
netically doped topological insulators [5T] . The monopole 
charge of Weyl nodes in these systems is ±1. Neverthe¬ 
less, HgCr 2 Se 4 [ 22 ], SrSi 2 [23] are expected to host Weyl 
nodes with monopole charge ±2, dubbed as the double- 
WSM. The topological invariant and enclosed Berry flux 
in double-WSM is twice that in a WSM, and conse¬ 
quently the one-dimensional chiral surface states in the 
former system possess a two-fold degeneracy. Although 
electro-magnetic responses in Weyl materials are rea¬ 
sonably well understood [6], stability of incompressible 
topological semimetals in the presence of quenched ran¬ 
domness is yet to be explored and settled. This is the 
quest that has recently culminated in a surge of analyti¬ 
cal [U [211 I26U34] and numerical [3] [35l [371444] works and 
we pursue it here for WSM and double-WSM, using nu¬ 


merical and analytical methods. 



FIG. 1: (Color online) MDOS in (a) WSM (N w = 1) and 
(b) double-WSM. WSM remains stable up to W c = 3.3 ±0.1, 
beyond which the mean DOS at e = 0, f)(0) is finite, and the 
system becomes a CDM. The double-WSM visibly turns into 
a CDM for weak enough disorder W > 1.0. 

We address the stability of these two systems against 
random quenched charge impurities by analyzing the 
mean density of states (MDOS) at the zero energy, where 
non-degenerate valence and conduction band touch each 
other. Our central results are: (i) For sufficiently 
weak disorder, while WSMs remain stable, double-WSM 
undergoes a BCS-like weak coupling instability toward 
the formation of a compressible diffusive metal (CDM), 
where the MDOS at zero energy is finite [see Fig. [l]. 
(ii) WSMs undergo a disorder driven quantum phase 
transition (QPT), beyond which the system becomes a 
CDM. (iii) Across the WSM-CDM transition MDOS dis¬ 
play single-parameter scaling and within our numerical 
accuracy the critical exponents at such itinerant quan¬ 
tum critical point (QCP) appear to be insensitive to the 
number of Weyl pairs (Nw) (see Table |T]). 

Model: A paradigmatic two-band toy model 

Hw = [iVi(k)ai + 7V 2 (k)(j2 + 7V 3 (k)cr 3 ] (1) 

k 

can describe different members of the Weyl family, where 





FIG. 2: (Color online) Finite energy phase diagram of a dirty 
(a) WSM (N w = 1) and (b) double-WSM, for L = 220. 
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<y are standard Pauli matrices. Fermionic annihila¬ 
tion operators c s ,k with spin-projections s =f, 3 and 
wave vector k, constitute the two component spinor 
= ( c t,k’ c t,k)- A WSM is found upon choosing 
Nj( k) =tsm(kja) for j = 1, 2, N 3 ( k) = N%( k) + IVf (k), 
where 2 iV 3 (k) = tcos(k 3 a) and a is the lattice spac¬ 
ing. In this model, the number of Weyl pairs (Nw) 
can be tuned efficiently by the Wilson mass 7V|(k) = 
t'[b — cos(fcia) — cos(/c 2 a)]- A double-WSM can be 
constructed by taking IVi(k) = ti[sin(/cia) — sin(fc 2 a)], 
N 2 (k) = ti cos(fcia) cos(fc 2 a), and -Af(k) = t '[2 — 
sin(fc x a) — sin (k y a)], while keeping (k) unaltered. We 
implement these tight-binding models on a cubic lattice 
with periodic boundary in each direction pEBj . 

Disorder: The quintessential properties of dirty Weyl 
systems can be established from their effective low energy 
theory in the close vicinity of the Weyl points. The low 
energy Hamiltonians for WSM and double-WSM are 


ffi 

H 2 


Tj: [-iv (cr i d\ + a 2 d 2 + ra 3 d 3 ) + V(r)] T t , (2) 
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respectively, where v ~ ta, mT 1 ~ t\a 2 , r = ± represent 
left and right chiral sectors respectively, and we set t = 
t' = t\ = 1 = a. Effect of random impurities is captured 
by V(r), distributed uniformly and independently within 
[— T’ T"]> and the MDOS is numerically evaluated using 
the kernel polynomial method nano. 

To gain insight into the role of disorder in these sys¬ 
tems, we can perform disorder averaging, assuming a 
Gaussian white noise distribution with zero mean, i.e., 
{{V (r)F(r'))) = A<5 3 (r — r') and arrive at the replicated 
Euclidean action 

S n = J d 3 xdt 

- d^dtdt'^a)^^^)^, ( 4 ) 


where a, b are replica indices and H n corresponds to the 
Hamiltonian from Eqs. |2| and ([3J) in the clean limit. The 
scale invariance of physical observables (v and m) dic¬ 
tates the following space-time(imaginary) scaling ansatz: 
(x, y) — > e l / n (x, y), z —► e l z and t —I e l t, accompanied by 


FIG. 3: (Color online) MDOS at zero energy £>(0) in double- 
WSM as a function of l/W. Inset: scaling of g(0) with W. 


the rescaling of fermionic field 4/ —> e (™ + 2 ) ( \jq where 
l ~ log F is the scaling parameter. The scaling dimen¬ 
sion of disorder coupling is [A] = 1 — A _ Hence, suffi¬ 
ciently weak disorder is an irrelevant (since [A] = — 1) 
and a marginally relevant (since [A] = 0) perturbation 
in WSM and double-WSM, respectively [351. Therefore, 
WSM (double-WSM) is expected to be stable (unstable) 
in the presence of sufficiently weak randomness. 

Notice that [A] = 2z — d [31], where d is dimen¬ 
sionality of the system and 0 is the dynamic critical 
exponent, together governing the scaling of mean DOS 
g(e) ~ |e| d / z_1 . Therefore, with z = 1 and §, g(e) ~ |e| 2 
and |er|, respectively for WSM and double-WSM, in agree¬ 
ment with our numerical findings, see Fig. [Tj 

Stability: WSM evidently remains stable for weak dis¬ 
order (W < 3.3) and MDOS at zero energy g( 0) = 0, in 
agreement with our scaling theory [see Fig. [lja)]. How¬ 
ever, for strong disorder WSM appears to undergo a QPT 
and enter into the CDM phase, where g( 0) becomes finite. 
A finite energy phase diagram of a dirty WSM is shown 
in Fig. [2j[a) . This observation is in qualitative agreement 
with various field theoretic [Tj [211 f21jf25] , and numerical 
analyses for three dimensional Dirac [31 [35U381 SO] and 
Weyl [351 HQ J2: semimetals. We will discuss the nature 
of such QPT in a moment. 

The scaling analysis suggests a BCS-like instability of 
double-WSM toward the formation of CDM for infinites¬ 
imal randomness, and a phase diagram of this system 
is shown in Fig. §b). By contrast, available data of 
p(0) for double-WSM in a finite system suggest a pu¬ 
tative threshold value of disorder (Wth ~ 1), only be¬ 
yond which p(0) is visibly finite [Fig. |TJ[b)] . To exam¬ 
ine whether the observed fmite-g(O) is a consequence of 
a W = 0 + instability, we compare g( 0) vs. l/W (see 
Fig-i- It is evident that larger systems are required 
to pin the exponential onset of g( 0) for sufficiently weak 
W, as W t h ~ l/log(A). Still, g(0) depicts a overall good 
agreement with exponential decrease with l/W. In addi¬ 
tion, p(0) fits very well with the celebrated BCS scaling 
form p(0) ~ exp(—A/W), with non-universal parameter 
A = 7.2 ± 0.8 [see the inset of Fig. [3], Instability of 
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N w 

W c 

z 

vm 

Vw 

v L 

1 

3.3(0.1) 

1.42(0.05) 

0.97(0.1) 

0.72(0.2) 

0.95(0.1) 

2 

2.5(0.1) 

1.38(0.05) 

1.1(0.15) 

0.72(0.2) 

1.1(0.15) 

4 

2.2(0.1) 

1.49(0.05) 

0.86(0.06) 

0.8(0.15) 

0.9(0.1) 


TABLE I: Comparison of critical disorder for WSM-CDM 
QPT (Wc), dynamic critical exponent ( 2 ), and correlation 
length exponent (is) extracted from the scaling of MDOS 
(see text) for WSM with Nw = 1,2,4 [45] • Quantities in 
parentheses denote the fitting error. Near WSM-CDM QPT 
at strong disorder all nodes get coupled and with increasing 
back-scattering channels or Nw, W c gradually decreases. 


double-WSM against weak enough disorder is, however, 
insensitive to the nature of disorder, and similar out¬ 
come holds for magnetic disorders 05]. But, g(0) starts 
to deviate from BSC scaling for W > 5.0 and falls below 
the exponential line. Such behavior can be attributed 
to the well-known Anderson transition of a three dimen¬ 
sional metal, across which MDOS does not display criti¬ 
cal behavior, but decreases monotonically [38]. We antic¬ 
ipate that for W > 5.0 the double-WSM falls within the 
basin of attraction of metal-Anderson insulator critical 
point. However, due to a logarithmic onset of a metallic 
phase in dirty double-WSM, for sufficiently weak disor¬ 
der and/or small system size, quasiparticle excitation can 
retain their ballistic nature over a large energy scale [see 
Fig- [Hb)]. 

Criticality: Now we investigate the scaling of MDOS 
across disorder driven QPT in WSM. Total number of 
states AT(e, L) below the energy e in a system of linear 
size L is proportional to L d , and in general is a function 
of two dimensionless variables L/£ and e/sq. While the 
correlation length diverges as £ ~ 8~ v , the correspond¬ 
ing energy scale vanishes according to £q ~ 8 UZ , as one 
approaches the QCP (S —> 0), where 8 = (W — W c ) jW c 
measures the deviation from the QCP ( W c ) and v is the 
correlation length exponent 03 08] . Consequently, 

M(e,L) = (L/Z) d g(e6- v *,L 1 /''5), ( 5 ) 

where Q is an unknown scaling function. From the defini¬ 
tion of MDOS g(e, L) = L~ d cU\f (e, L)/de, we then arrive 
at the following scaling ansatz 

e {e,L)=8^ d - z ^(\e\5-" z ,L^5), ( 6 ) 

after accounting the particle-hole symmetry, g(e, L ) = 
g(—£, L), where T is also an unknown, but universal scal¬ 
ing function. Below we demonstrate the scaling analysis 
of MDOS in WSM with Nw = 1 and the results for 
Nw = 1,2 and 4 are summarized in Table [T] 05] . 

First we consider a sufficiently large system (L = 
220), so that the finite size effects are negligible and In¬ 
dependence in Eq. ([6]) can be ignored. From the scaling 
of MDOS at zero energy g(0) with disorder we estimate 
critical disorder for WSM-CDM QPT W c = 3.3 ± 0.1 
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FIG. 4: (Color online) (a) MDOS at zero energy £>(0) vs. 
disorder; (b) MDOS g(s) vs. e for W = 3.2,3.3, 3.4; (c) g(0) 
vs. 8, where 5 = (W — W c )/W c ", (d) c(8)~ 1 vs. <S, with c(<5) = 
for WSM with Nw = 1, assuming W c = 3.3 [45 .. 

[Fig.ga)]. At the QCP (8 = 0), the ^-dependence of g(e) 
must cancel out, demanding F(x) ~ a:" 1 , and therefore 
g(e) ~ |e|^ _1 . From Fig. gb) we obtain z = 1.42 ±0.05. 

In the metallic phase, the MDOS at zero energy g(0) is 
finite and serves as the order parameter. In this regime 
g(0) ~ S^ d ~ z ^ L/ and one can identify (d— z)v as the order 
parameter exponent (/3). However, such power law de¬ 
pendence of g(0) is valid when (<L. Therefore, we fit 
p(0) as 8 < ' d ~ z ' >1 ' for 8 > 0.06, and obtain vm = 0.97 ± 0.1, 
where vm is the correlation length exponent extracted 
from the metallic phase [see Fig. S')]- 

In WSM, the mean DOS scales as g(e) ~ c(S) 1 |e| d 1 , 
so that we recover g(e) ~ |e| 2 for d = 3, where c(8) ~ 
£(z-i)dv w an( j Vw j g t j le correlation length exponent ex¬ 
tracted from the WSM phase. However, it should be 
noted that for W < W c , mean DOS displays a smooth 
crossover from |e| 2 (for small e) to |e| (for large s) depen¬ 
dence. Therefore, estimation of vw depends crucially 
on the range over which we attempt to fit g(e) ~ |e| 2 , 
and accuracy of vw can be questioned. Nevertheless, 
by fitting the coefficient of |e| 2 with c(8)~ 1 , we obtain 
v w = 0.72 ± 0.2 [see Fig. gd)]. 

Data collapse: We now demosntrate that g(e) display a 
single parameter scaling across WSM-CDM QPT. First, 
we compare g(£)S~^ d ~ z ^ l/ vs. \£\8~ zv for L = 220. Ne¬ 
glecting the high energy part of spectrum (|e| > 0.5) out¬ 
side the weyl cones and extreme small energy (|e| < 10~ 2 ) 
where numerical accuracy is small, we find that all data 
from Fig. ga) collapse onto two separate branches, asso¬ 
ciated with the CDM and WSM phases [see Fig.ga)]. 

Next we delve into the finite size data collapse for 
MDOS at e = 0 and estimate u independently. Setting 
e = 0 in Eq. (g), we obtain g(0,L) = L z ~ d J : \t),8L 1 / v ). 
An excellent data collapse is achieved by comparing 
g(L, 0 )L d ~ z with 5L 1 /" for several systems with 80 < 
L < 180, W c = 3.3 and z = 1.42 [see Fig. gb)]. The cor- 

















FIG. 5: (Color online) (a) Single parameter scaling of MDOS 
in WSM (Nw = 1) with L = 220. Top (bottom) branch 
corresponds to CDM (WSM). (b) Data collapse for g(0,L) 
when Nw = 1. For large <5 deviation from single-parameter 
scaling stem from the Anderson transition at strong disorder. 


relation length exponent extracted from the best quality 
data collapse is vl = 0.95 ± 0.1. Thus, p(0) displays a 
single-parameter scaling and serves as an bonafide order 
parameter across the WSM-CDM QPT. 

Critical regime'. The crossover from the quadratic (at 
small energy) to the linear (for higher energies) scaling 
of p(e) allows us to estimate the crossover boundary be¬ 
tween the WSM and critical regime at finite energy when 
W < W c . For sufficiently weak disorder, g(s) ~ \e\ 2 over 
a wide range of energy. As the randomness is gradually 
increased (but still W < W c ), more and more degrees 
of freedom need to be integrated out (in the spirit of 
renormalization group) to wash out the effect of disor¬ 
der from the system. Consequently, the energy window 
over which g(e) ~ |e| 2 gets reduced and the region where 
g(e) ~ |er| increases, with increasing disorder. With this 
notion we numerically estimate the crossover boundary 
between the WSM and the critical regime at finite en¬ 
ergy [see Fig. [2ja)] . When W = W c , the mean DOS 
displays a |e|-linear dependence over the entire energy 
range (|e| < 0.5). For W > W c , |e|-linear behavior of 
g(e) ceases at finite energy, defining the boundary be¬ 
tween the CDM and critical regime [see Fig. |2ja)] . 

Conclusions: To conclude we show that WSM is sta¬ 
ble against weak disorder, but undergoes a QPT and 
becomes a CDM at strong disorder. Across such QPT 
MDOS display single-parameter scaling. The critical ex¬ 
ponents (u, z) appear to be independent of the num¬ 
ber of Weyl nodes (Nw) [see Table [I], since contribu¬ 


tion from any fermionic-bubble vanishes in the vanishing 
replica limit [28]. The extend of the critical regime at fi¬ 
nite energy associated with such QCP [see Fig. |2ja)] can 
also be measured from ARPES or scaling of specific heat 
(C v ~ T d / Z ) in various WSMs CM] and topological 
Dirac semimetals (two superimposed copies of WSMs), 
such as Cd 2 As 3 JJS, Na 3 Bi m- In contrast, the double- 
WSM becomes CDM for weak (infinitesimally small in 
the thermodynamic limit) disorder. While in the metal¬ 
lic phase C v ~ T, in WSM and double-WSM specific 
heat scales as C v ~ T 3 and T 2 , respectively [J5]. There¬ 
fore, as a function of temperature specific heat in double- 
WSM should display a smooth crossover from T 2 to In¬ 
dependence as temperature is gradually decreased [see 
Fig- §b)]. Generalization of scaling analysis dictates 
that weak disorder is a relevant perturbation in triple- 
WSM (monopole charges ±3), since [A] = 1 — ^ | 

for n = 3. Therefore, among various three-dimensional 
topological semimetals only conventional WSM is stable 
against weak disorder m- 

Finally, we discuss the transport phenomena in disor¬ 
dered Weyl systems. In weakly disordered WSM (W < 
W c ) the optical conductivity (in collisionless regime) dis¬ 
plays a smooth crossover from a n (Q) ~ Q to Q 1 / 2 de¬ 
pendence as frequency (Q) is increased, closely following 
the phase diagram in Fig. [2ja) for j = x,y,z [34] • In 
strong disorder regime ( W > W c ) <Jjj (Q) becomes finite 
as fl —> 0. By contrast, in double- and triple WSMs 
cr zz (C.) ~ H, while cr xx / yy (Q) ~ Q 1 / 71 at high frequency 
with n = 2 and 3, respectively. However, as Cl —> 0, <jjj 
becomes finite in these two systems for arbitrary strength 
of disorder and for any j. Scaling of dc conductivity 
(collision dominated) follows the ones those for optical 
conductivity upon taking Q —>• T. Thus, in future one 
can probe the transport properties to establish the global 
phase diagram of disordered Weyl materials at finite fre¬ 
quency and temperature. 
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The Supplementary Materials contain: 

1. The renormalization group analysis and additional numerical results for disordered double-WSM, 

2. Computation of scaling dimension of disorder coupling for arbitrary quasi-particle dispersion in a d-dimensional 
Weyl system, 

3. The numerical analysis for WSM with Nw = 2 and 4, where Nw is the number of Weyl pairs, 

4. Description of tight-binding model for WSM and double-WSM, 

5. Some essential details of numeric methods, which we employed to analyze the behavior of MDOS in dirty WSM 
as well as double-WSM. 

6 . Scaling of specific heat (C v ) with temperature (T) for a dirty double-WSM. 

• RG analysis and additional numerical results for double-WSM: While scaling analysis suggest that 
disorder is a marginal perturbation in double-WSM, to establish that weak disorder is marginally relevant one needs 
to account for the quantum corrections. The replicated action (imaginary time) in the presence of generic disorder 
reads as 

/ a 2 _ n2 o A f) A n f 

d 3 xdt [dt + 0-1 2 2m 1 - 0 - 2 -^^ - ivTa 3 d z \^a(x,t) - -y j dAxdtdt' (*!*&)(*,t') 

3 

- f d 3 xdtdt' (Vi<7iV a ) {Xtf) {^l<Ti^ b ) {x>t>) , (7) 

3 =1 J 


where a, b are replica indices. Here we allowed all types of disorder with distinct bare strengths to examine the low 
energy behavior in dirty double-WSM in the presence of disorder of arbitrary nature. 


Upon integrating out the fast Fourier modes with Ae 1 < + v 2 k1 < A, where A is the ultraviolet cut-off for 

energy, we arrive at the following renormalization group (RG) flow equations to quadratic order in disorder coupling 


= Aq + -Ao (Ai + A 2 ) + A 1 A 2 + A 0 A 3 , = 9 [Ag + A 2 + 2 A 1 A 2 + 3 A 1 A 3 + A 2 + 2 A 0 A 2 — AqAi] , 

= 2 [^0 + A 2 + 2 A 1 A 2 + 3 A 2 A 3 + A 2 + 2AqAi — A 0 A 2 ] , = (Ai + A 2 ), ( 8 ) 

after taking Ajm/(2nv) —> A j for j = 0,1, 2, 3. Notice that disorder coupling A 3 does not get generated through coarse 
graining if the bare model is deviod of such elastic scatterer. Therefore, a bare model with A 3 = 0, remains closed 
under RG. However, both Ai and A 2 gets generated from A 0 , thus these three couplings needs to be simultaneously 
accounted for to keep the model closed under RG. Interestingly if we start with a bare model with only Ao 7 ^ 0, 
A i(l) = A 2 ( 1 ) for all RG time (/), and one arrives at a simpler set of flow equations 


dAo 

~df 

dA 2 

dl 


— Ag + A 2 ± + 3A 0 Aj_, — 2 (Ag + 3A^ + A 0 A_l) , (9) 

at the one-loop level, where Ai = A 2 = Aj_. The above set of flow equations [see Eq. (| 8 |] supports only one unstable 
fixed point at (A 0 , A 1; A 2 , A 3 ) = (0,0,0, 0). Thus weak disorder of arbitrary nature is always a marginally relevant 
perturbation in double-WSM and drives the system immediately into the CDM phase. 


To test the validity of our RG analysis we performed numerical analysis of MDOS at zero energy for all four 
different types of disorder. Result for regular potential disorder (A 0 ) is shown in the main part of the paper. In this 
Supplementary Material we present the numerical results for MDOS at zero energy in double-WSM in the presence 












FIG. 6: (Color online) (a) Mean density of states (MDOS) at zero energy in double-Weyl semimetal in the presence of various 
types of magnetic impurities, denoted by the Pauli matrices [see also Eq. 01- Here, IF denotes the strength of disorder. Notice 
that MDOS at zero energy falls on the same exponential curve irrespective of the nature of impurity scatterer. Here, we used 
the fitting function g( 0) = ao + ai exp (— 02 IF -1 ), with ao = 8.6 x 10 -4 ± 9.06 x 10 -5 , a 1 = 3.8 ± 0.76, 02 = 1.64 ± 0.59. (b) 
BCS scaling of MDOS at zero energy £>(0) vs. IF, g(0) ~ exp (—A/IF), with A = 1.65 ± 0.597. 



FIG. 7: (Color online) Flow of disorder couplings for (a) Ao(0) = 0.1, and Ai(0) = A2(0) = A 3 (0) = 0; (b) Ai(0) = 0.05 
and A o (0) = A 2 (0) = A 3 (0) = 0; (c) A 2 (0) = 0.05 and A o (0) = Ai(0) = A 3 (0) = 0; (d) A 3 (0) = 0.05 and A o (0) = Ai(0) = 
A 2 (0) = 0, where Aj(0) represents the bare value of dimensionless disorder coupling (see text). Here, l is the RG time. 




























of three different types of magnetic disorder Ai, A 2 and A 3 , and results are shown in Fig. [bj The fact that ADOS 
at zero energy for arbitrary magnetic disorder (A j for j = 1,2,3) falls on one curve is quite interesting and can be 
substantiated from the solution of the RG flow equations, displayed in Eq. Q. As shown in Fig. [7j irrespective of the 
initial conditions, a particular disorder coupling Aq (associated with potential disorder/random charge impurities) 
always diverges first, leading to BCS-like instability of double-WSM for arbitrarily weak disorder of any nature. 
Consequently, the MDOS at zero energy for all magnetic disorder falls on the same curve as shown in Fig. [ 6 j It is 
worth mentioning that irrespective of disorder coupling, we always find that A i(l) = A 2 (Z) for l 1, as shown in Fig. [7] 

To test the rotational symmetry of double-WSM in the x — y-plane in the presence of any disorder, we next compute 
the two-point correlation function or the self-energy correction to the leading order, which is given by 

k) = (— iu>) (A 0 + Ai + A 2 + A 3 ) - (10) 

that only depends on external frequency (w). Thus, rotational symmetry of double-WSM in the x — y plane remains 
unaffected in the presence of arbitrary disorder, and the system becomes unstable toward a CDM phase for arbitrarily 
weak disorder. 

Finally, to establish the BCS-scaling of MDOS one needs to find the self-consistent solution of the quasi-particle 
scattering life-time (r) at zero energy, since MDOS at zero energy p(0) is inversely proportional to r [Tj. By employing 
self-consistent Born approximation we arrive at |2] 

= 1 P (“S)’ (U) 

after substituting g(e) ~ |e|, where A is a non-universal model and material dependent constant, and A is the 
high-energy cut-off. Thus g(0) displays BCS-scaling as we found in our numerical analysis and predicted from the 
RG analysis. 

Notice that here we considered only the intra-valley disorder potentials. For Anderson transition one needs to 
account for the inter-valley scattering as well. In a numerical analysis such inter-valley scattering is always present. 
However, inclusion of inter-valley scattering does not alter the instability of ballistic fermions in double-WSM toward 
the formation of a metallic phase for arbitrary weak disorder. 

• Scaling dimension of disorder: We now demonstrate the scaling dimension of disorder coupling in a d- 
dimensional system with arbitrary quasi-particle dispersion, captured by the low-energy Hamiltonian H( k). The 
imaginary-time action for such system after performing the disorder averaging is given by 

S= j d 3 xdt (*l[d t + H(k^ -*V)]*a) (Xit) - J J d 3 xdtdt'(^ a )^ t) (^ b ) {x ^, (12) 

where a and b are replica indices. We here perform a slightly different rescaling of space-time(imaginary) coordinates 
from the one shown in the main part of the paper, under which (x, t ) —» (e ; x, e zl t) and the parameter z is chosen such 
that all physical parameters in H(k —» —iV), such as Fermi velocity (v), quasi-particle mass (to), remain scale invari¬ 
ant. The parameter z also determines the scaling of MDOS according to g(e) ~ \s\ d / z ~ 1 . The replicated action S then 
remains scale invariant when the fermionic field is rescaled according to T —» e ~ dl / 2 ^. Under such rescaling of coordi¬ 
nates and field, the disorder coupling A —> e( d ~ 2z ^ l A. Therefore, scaling dimension of disorder coupling is [A] = 2 z — d. 

From the main part of the paper, where we performed an anisotropic scaling over various spatial coordinates 
(to keep v and m scale invariant) with t —> e l t, we know that [A] = —1 for WSM and 0 for double-WSM. 
Therefore, for these two systems z = 1 and 3/2, respectively. As a result, the MDOS in these two systems 
respectively scales as g(e) ~ |er| 2 and |e|. These results are announced in the main part of the paper and also con- 
sitent with numerical computation of MDOS in WSM and double-WSM, shown in Fig. 1 of the main part of the paper. 

• Numerical analysis for WSM with N\y = 2 and 4; We now present the numerical analysis for W c (critical 
disorder for WSM-CDM transition), and computation of various critical exponents (v and z) in WSM with Nw = 2 
and 4, where Nw is the number of Weyl pairs. The results are already quoted in Table I of main part of the 
paper. These quantities are extracted by using identical methods, described in the main paper for Nw = 1. The 
corresponding analysis for Nw = 2 and 4 are presented in Fig. [ 8 j The data collapse for MDOS at finite energies in 
WSM and CDM phases, as well as that for MDOS at zero energy in the CDM phase for Nw = 2 and 4 are shown in 




FIG. 8: (Color online) (a) Mean DOS at zero energy £>(0) vs. disorder (IF), mean DOS g(s) vs. e for W = 2.4, 2.5, 2.6, £>(0) 
vs. 8, where $ = (W — W c )/W c , and c(5)“ 1 vs. (5, where c(<5) = , for WSM with Nw = 2 (from left to right, and top 

to bottom), (b) Similar quantities for Nw = 4. Extracted values of W c , z and v are quoted in the figures. 



FIG. 9: (a) Single-parameter scaling of MDOS at finite (left) and zero (right) energy in WSM with Nw = 2. (b) Similar 

quantities for Nw = 4. Top (bottom) branch in the left panel of (a) and (b) corresponds to metal (Weyl semimetal). 


Fig. [9] 


• Tight-binding model for WSM and double-WSM: The tight-binding model that support a WSM, is 
described by the following momentum dependent form factors, 7Vi(k) = tsin(fcia), A^k) = tsin(fc 2 a), N£( k) = 
tcos{k 3 a)/2 and 7V|( k) = t'[b — cos(kia) — cos(fc 2 a)], where a is the lattice spacing and N 3 (k) = N 3 (k) +7Vf(k). The 
Hamiltonian is then given by 

H w = J2^l [JVi(k)tri + N 2 (k)a 2 + N 3 (k)a 3 ] tf k , (13) 

k 


which on a cubic lattice with periodic boundary in each direction translates into the follwoing tight-binding model 


* = £ 


- ’Fjaa’J'r+es + £ *1 

it t! 

2 a > ~ 2 aa 

^r+£j + H.C 

L J=l,2 




£*J [2tV 3 + V(r)) * P , 


(14) 


where ’Fj = (c ri |, c ri j_), and c FiS is electronic annihilation operator at site r with spin projection s =t,i- Nearest- 
neighbor sites are connected by the unit vectors Cj, for j = 1,2,3. For b = 2, a pair of Weyl nodes are located at 
k = (0, 0, ±1)2^, while for 6 = 0, two pairs of Weyl nodes are found at k = (1, —1, ±1)^ and (—1,1, ±l)g= Finally, 
when t' = 0, four pairs of Weyl nodes are found at k = (0,0, ±1)^, (1, — 1, ±1)^, (—1,1, ±1)^ and (—1, —1,±1)|= 
Therefore, by tuning various parameters of the above tight-binding model one can realize WSMs with different Nw 
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FIG. 10: (Color online) Scaling of C v /T with temperature T (measured in units of t ) in dirty double-WSM. Strength of disorder 
is quoted in the figure. 
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As announced in the main part of the paper that by setting Aq(k) = tj^sinffeia) — sin(fc 2 a)], N 2 (k) = 
ti cos(fcia) cos(£; 2 a), 7Vg(k) = tcos(k 3 a ) and N 3 (k) = t'[ 2 — sin(feia) — sin(fc 2 a)], one can realize a double-WSM 
in the vicinity of k = i^). In cubic lattice with periodic boundary in each direction, the corresponding 

tight-binding model is 


H 


E 


r+ e 3 + | 


E*J {2t'a z + V(r)\ * r . 


^ai^r+e,] 


7 E K+«e 2 4E + H.C 

a=± j=1,2 J 

(15) 


Here, V (r) captures the effect of random quenched charge impurities. 


• Numerical method: We address the stability of various members of the Weyl family by analyzing the MDOS at 
zero energy. The MDOS is numerically calculated using the kernel polynomial method [SI [!]• The MDOS is evaluated 
using the following definition 

g(e) = Tr S(e - H), (16) 

where Tr is evaluated stochastically and the delta-function is exapnded using Chebyshev polynomial. We usually 
take 4096 Chebyshev moments and few (~ 8) trace vectors to calculate MDOS. We ensure the convergence of MDOS 
with respect to number of Chebyshev moments and number of random vectors. We typically take ~ 30 disorder 
realization to further minimize the residual statistical error, which is quite small in large system (L = 220) due to 
the self-averaging nature of MDOS. 


• Specific heat (C v ) in dirty double-WSM: The specific heat is related to the MDOS according to 


c„ 


r de g{£) g2 

J- oo 4 cosh“(/3e/2) 


(17) 


where /3 = 1 /(fcgT) is the inverse temperature and we set ks- In numerical calculation of C v from the tight-bonding 
model, integral over energy is restricted by the band width. 


Since in clean double-WSM MDOS g(e) ~ |e|, the specific heat C v ~ T 2 . On the other hand, in the metallic phase 
C v ~ T. To detect the onset of an metallic phase, we thus compare C v /T vs. T for a wide range of temperature. The 

For any strength of disorder (W) C v /T becomes finite as T —>• 0. For weak enough 
T for a wide range of temperature, while such range of T gradually decreases with increasing 


results are displayed in Fig. flO 
disorder C v /T 


disorder in the system. Therefore, in the presence of even infinitesimal disorder a double-WSM turns into a diffusive 
metal at lowest energy scale. However, over a wide temperature regime quasi-particles retains their ballistic nature 
as reflected from C v /T ~ T behavior, which, however, diminishes with increasing strength of disorder. 




























From Fig. 10 we also note that with increasing strength of disorder the magnitude of specific heat increases monoton- 
ically, stemming from the fact that the MDOS in the system with increasing randomness. Therefore, with increasing 
disorder system becomes more metallic. 


[1] R. Shindou, and S. Murakami, Phys. Rev. B 79, 045321 (2009). 

[2] P. Goswami, A. H. Nevidomskyy, Phys. Rev. B 92, 214504 (2015). 

[3] K. Kobayashi, T. Ohtsuki, K-I. Imura, I. F. Herbut, Phys. Rev. Lett. 112, 016402 (2014). 

[4] A. Weifie, G. Wellein, A. Alverman, and H. Feshke, Rev. Mod. Phys. 78, 275 (2006). 




